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Abstract. We present an illustration of using a quantum three-body code being prepared for public release. The code is based 
on iterative solving of the three-dimensional Faddeev equations. The code is easy to use and allows users to perform highly- 
accurate calculations of quantum three-body systems. The previously known results for He3 ground state are well reproduced 
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INTRODUCTION 

The quantum few-body problem is important for investigating physical processes at practically all possible length and 
energy scales. For instance, three-body models can be employed for describing nuclear reactions [1, 2, 3, 4], electron- 
and positron-atom collisions [5, 6], and chemical reactions [7]. The developments of the last decade demonstrated 
the importance of three-body processes for understanding the dynamics of ultra-cold gases [8]. An ability to solve a 
few-body problem directly would also be beneficial to theorists for testing, for instance, effective-field theories [9, 10]. 

The three-body problem, however, has sufficient intrinsic complexity that it often inhibits or prevents non-experts 
in few-body calculations from considering realistic three-body models and from employing physically correct repre- 
sentations [11]. Accordingly, a standard, easily operable and rigorously constructed tool for three-body calculations 
will be beneficial for a broad physical community. Such a tool should be tested independently to ensure its usability 
and applicability. This work is a result of a collaboration between the authors of this tool, developed at the University 
of Kentucky, and a research group in JINR performing independent tests. 

In the following sections we describe the equations being solved and report the results of the tests we have 
performed. 

FORMALISM 

The three-body code being tested is based on solving Faddeev equations in configuration space. The complete and 
mathematically rigorous theory of Faddeev equations can be found in books on the topic [12, 13, 14]. Here we only 
sketch out the gross features important for understanding and using the code. 
We start by describing the physical model from the three-body Hamiltonian 

H=H + V 3b {x i ,y i )+'£V i (x i ), (1) 

i 

where Ho stands for the kinetic energy of the three particles, V, (x,) is the interaction potential acting in the pair ;', 
and V3b(x,,y;) is a short-range three-body interaction. (In the following description the latter will be omitted only 
for simplicity. Taking into account the three-body interactions, however, does not produce any practical or principal 



difficulties.) The configuration space of the three particles is described in terms of 3 sets of Jacobi coordinates 
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The set of coordinates ;' describes a partitioning of the three particles into a pair (jk) and a separate particle i. Faddeev 
decomposition represents the wave function *P in terms of a sum over all possible partitioning of the three-body system 

V = £*,(x*,y i )- (3) 

i 

Faddeev components, 4>,, satisfy the following set of equations [12] 

(-A x -A y + Vi-(x ( -) -E)4>,(x,,ji) =VJ(x ( )£*t(x t ,yt) , (4) 

where x, and y, are mass-weighted Jacobi coordinates, Vj is the interaction potential in the ;-th pair and E is the total 
energy of the system. It is not difficult to prove that the exact wave function of the three body system can be uniquely 
constructed from the Faddeev components by means of Eq. (3). 

The equations in six-dimensional space can hardly be solved directly, and some partial analysis is necessary. We 
consider the states with zero total angular momentum. The angular degrees of freedom corresponding to collective 
rotation of the three-body system can be separated [15] and the kinetic energy operator reduces to 

where x, y and z are so called intrinsic coordinates 

jc=|x|, y = |y|, z=^, x,ye[0,oo), z e[-l,l] . (6) 
xy 

In the case of identical bosons Faddeev components take identical functional form, which makes it possible to reduce 
the system of three equations (4) to one equation 

(H + V(x)-E)(j>(x,y,z) = -V(x)P(l>(x,y,z), (7) 

where 
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Assuming that in each two-body subsystem only one bound state exists, we can write the asymptotic boundary 
conditions for the Faddeev component as follows 



<j)(x,y,z) ~ 9i(x)e ^ +A{x/y,z) {x 2 +y2 y/A ' (9) 



where (p2(x) stands for the wave function of the two-body subsystem bound state, k y = \fE~2—E~i, = ^J— £3, E2 is 
the two-body bound state energy, and £3 the energy of the three-body system. For three-body bound states the first 
term corresponds to virtual decay into a particle and a two-body bound system, while the second term corresponds to a 
virtual decay with an amplitude A (x/y,z) into three single particles. The term corresponding to the latter configuration 
can generally be neglected for the states below the three-body threshold. Therefore, at sufficiently large distances R x 
and R y , the asymptotic boundary conditions for the Faddeev component are 
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(10) 



For bound state calculations Dirichlet or Neumann boundary conditions can also be employed. 

The important property of the Faddeev components which makes them suitable for numerical solution is their simple 
asymptotic form. For instance, each of the Faddeev components holds only bound states of the corresponding two- 
body subsystem. In this respect the x, coordinate is the internal coordinate of the corresponding two-body cluster and 
the y, coordinate plays the role of a reaction coordinate for all the states below the 3-body (break-up) threshold. This 
simple physical meaning of the coordinates suggest a natural requirement for discretizing the corresponding degrees 
of freedom: the discrete analogs of the x, coordinate should reproduce the spectrum of the ;-th cluster correctly, and 
discrete analogs of the y, coordinate must describe the scattering states reasonably well. These necessary requirements 
are easy to check prior to performing actual calculations, and they also provide a solid ground for a reasonable degree 
of automation for choosing the parameters of the numerical scheme. Another advantage of the Faddeev equations is the 
asymptotic decoupling of the components. The right-hand side of the equation (4) is, roughly speaking, exponentially 
small if the third particle is at larger distance than the typical size of the two-body bound state. This means that 
at longer distances |y, | > y max the Faddeev components rapidly decouple, and calculations can be performed in the 
regions as small as the size of the largest two-body subsystem bound state. 

The advantages of the Faddeev approach can be exploited even further when dealing with short-range interactions; 
i.e. assuming that the potentials V, are zero (or negligible) outside the region |x,| < x max . To clarify this, consider 
the component 4>, in the asymptotic region |y,| — > 00 , x, e supV,, where it satisfies a Schrodinger equation with the 
corresponding channel Hamiltonian 

where Hq is the free three-body Hamiltonian. This property of <t>; suggests that, rather than calculate <t>, directly, we 
instead calculate T, (Eq. 11), which is better localized in configuration space 
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These T, are non-zero only for small x, = |x,| < x max and small y,- = |y, | < y max - Accordingly, they are localized to 
a region that can be much smaller than the typical size of a two-body bound state. This feature leads to substantial 
computational savings [16]. The T,- satisfy the following integral equations 
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(12) 



where R2j(E) are the resolvents of the corresponding channel Hamiltonians. Since supr, C supV,, the T, are more 
suitable for numerical approximation than the original Faddeev components. Furthermore, if the equation is being 
solved using an iterative technique, then no explicit representation for the integral operators R2j(E) is required. In this 
case we only need to calculate the action of the integral operator on the T;, which can be done with high computational 
efficiency by numerically solving the corresponding differential equation with appropriate boundary conditions. We 
call this computational scheme a Localized Component Method (LCM). 



COMPUTER CODE 

In order to construct a discrete analogue of the system of equations (12) we employ quintic Hermite splines together 
with the orthogonal collocations method [17]. A detailed description of the procedure is given in [18]. This high-order 
method guarantees fast convergence with respect to the number of grid points, sparse matrix structure for the discrete 
analog of the equation (4), and fast calculation of matrix elements that makes it possible to avoid storing big matrices 
in computer memory. 
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FIGURE 1. A screenshot of the configurator. 



The code is written in Java and consists of two parts. The first part is a configurator that simplifies composing 
the necessary configuration files. The configurator allows the user to set masses of the interacting atoms, to specify 
identical particles in the system, to choose a potential model, to set the cutoff distances R x and R y and to set the 
number of grid points to be employed in the calculation. It also generates a mesh with L2-optimal point distribution 
which ensures the best possible approximation of the three-body wave function in the asymptotic region. In Fig. 1 we 
show an example of the configurator screenshot. 

The second part is the three-body computational kernel based on the LCM approach. The kernel is currently capable 
of three-body calculations below the three-body threshold with a limitation of no more than one two-body bound state 
contributing to each asymptotic channel. This includes bound states, elastic scattering and chemical reactions below 
the first vibrational excitation threshold. Typical computational time can take from minutes to hours, depending on the 
physical system and the size of the grid. 



RESULTS 

We apply the code to the calculation of binding energies of the Helium trimer 4 He3 three-atomic system, to verify 
that we can reproduce the known calculated properties of helium trimer ground and excited states. This system is 
very particular about the approach being used, as the trimer binding energy is extremely small and a large volume of 
the configuration space should be treated, but the interaction features very strong repulsion at short distances, which 
requires very precise numerical methods to be used. 

Experimentally, helium dimers have been observed for the first time in 1993 by the Minnesota group [19], and in 
1994 by Schollkopf and Toennies [20]. Later on, Grisenti et al. [21] measured a bond length of 52 ±4 A for 4 He2, 
which indicates that this dimer is the largest known diatomic molecular ground state. Based on this measurement they 
estimated a scattering length of 104^^ 8 A and a dimer energy of 1 . l^R | pi]- ^ n the latter investigation [22] the 
trimer pair distance is found to be 1.1 ^n'g nm in agreement with theoretical predictions for the ground state. 

Many theoretical calculations of these systems were performed for various interatomic potentials [23, 24]. Varia- 
tional [26, 25], hyperspherical [8, 27, 28] and Faddeev techniques [16, 18, 29, 30, 31, 32] have been employed in this 
context. It was found that the Helium trimer has two bound states of total angular momentum zero: a ground state of 
about 126 mK and an excited state of Efimov-type of about 2.28 mK. Experimentally this Efimov-type[33] excited 
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FIGURE 2. Convergence of the Helium trimer ground state energy on the grids of N x = N y points; N z = 5. 



TABLE 1. 
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state has not yet been observed (see, e.g., [34] and refs. therein). It should be mentioned, however, that the year 2006 is 
noticeable due the first convincing experimental evidence for the Efimov effect in an ultracold gas of Caesium atoms 
[35, 36]. 

In present calculations we employed the code based on the Faddeev differential equations (5) with bound- 
ary conditions (10). As He-He interaction we used the semi-empirical LM2M2 potential [23]. We use — 

4.0026032197 a.u.m for the mass of the 4 He atom and — ^-r 



7— = 12.1 1928 K A 2 , unlike many three-body Calcula- 
te 8 

tions, see, e.g., [34], where a rounded value of the coefficient has been used. 

Investigation of the bound state energy convergence with respect to the number of grid points demonstrates that 
even a moderate number of points in variables x and y is sufficient to get up to six accurate figures for the energy of 
the ground state (Fig. 2). 

The 4 He dimer binding energies, 4 He- 4 He scattering lengths and mean values of the radius < R > and V < R 2 > 
obtained with the LM2M2 potential [23] are shown in Table 1 in comparison with experimental data [21]. All the 
values agree with an experimental estimation of Ref.[21] within quoted errors. The scattering length 4c of the 
system is bigger than the range of the potential by an order of magnitude. All these features characterize the Helium 
dimer as the weakest, as well as the biggest, diatomic molecule found so far. Due to the fact that the energy of the 
dimer is so small, one should expect that the £4 He trimer indeed possesses the theoretically predicted state of the 
Efimov type (see, [33, 34]). 

In Table 2 the results of trimer binding energies calculations obtained with LM2M2 potential are summarized. The 



TABLE 2. Results for binding energies of the He3 trimer for LM2M2 potential. 
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binding energies of the 4 He trimer ground (£4 He3 ) and exited (E^^) states are presented. These results demonstrate 
good agreement between different methods and show that the code competes well even against variational methods. 
It should be mentioned that the energy estimates obtained with the code are non-variational, and further variational 
improvements of the results are possible. 

We are planning to continue testing the code within current applicability limits, including scattering calculations, 
systems of distinguishable particles and modeling clusters of other rare gas atoms. 
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